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In this paper, we study massive gravity in the presence of Born-Infeld nonlinear electrodynamics. 
First, we obtain metric function related to this gravity and investigate the geometry of the solutions 
and find that there is an essential singularity at the origin (r = 0). It will be shown that due to 
contribution of the massive part, the number, types and places of horizons may be changed. Next, 
we calculate the conserved and thermodynamic quantities and check the validation of the first law of 
thermodynamics. We also investigate thermal stability of these black holes in context of canonical 
ensemble. It will be shown that number, type and place of phase transition points are functions of 
different parameters which lead to dependency of stability conditions to these parameters. Also, it 
will be shown how the behavior of temperature is modified due to extension of massive gravity and 
strong nonlinearity parameter. Next, critical behavior of the system in extended phase space by 
considering cosmological constant as pressure is investigated. A study regarding neutral Einstein- 
massive gravity in context of extended phase space is done. Geometrical approach is employed to 
study the thermodynamical behavior of the system in context of heat capacity and extended phase 
space. It will be shown that GTs, heat capacity and extended phase space have consistent results. 
Finally, critical behavior of the system is investigated through use of another method. It will be 
pointed out that the results of this method is in agreement with other methods and follow the 
concepts of ordinary thermodynamics. 


I. INTRODUCTION 

Regarding experimental agreements of Einstein gravity (EN) in various area of astrophysics and cosmology, moti¬ 
vates one to consider it as an acceptable theory. In addition, adding a constant term A in the EN-Hilbert action may 
lead to agreement between the results of EN-A gravity with dark energy prediction. 

On the other hand, general relativity is consistent with interaction of massless spin 2 fields, in which related 
gravitons are massless particles with two degrees of freedom. Since the quantum theory of massless gravitons is 
non-renormalizable [l!], one may be motivated for modifying general relativity to massive gravity. In order to build up 
a massive theory with a massive spin 2 particle propagation, one can add a mass term to the EN-Hilbert action. This 
will result into graviton having a mass of m which in case of to 0, the effect of massive gravity will be vanished. 
A class of massive gravity theory in flat and curved background which leads to absence @ and existence [3j of ghost, 
have been investigated. Also, the quantum aspects of the massive gravity and a nonlinear class of massive gravity in 
ghost-free field 0] have been explored in Refs. SI- Generalization to nonlinearly charged massive black holes was 
done in Ref. [7:]. More details regarding the motivations of massive gravity is given in Ref. (§]. 

In this paper, we are interested in studying the nontrivial adS massive theory that was investigated in [Ici[. The 
motivation for this consideration is due to fact that graviton shows similar behavior as lattice in holographic conductor 
El- In other words, a Drude like behavior is observed for the case of massless graviton in this theory which makes 
the role of graviton similar to lattice. Another interesting subject for study in this theory is metal-insulator transition 
E|. Recently, charged massive black holes with consideration of this theory have been investigated in 0. The 
P — V criticality of these solutions and their geometrical thermodynamic aspects have been studied 0,0. Also, the 
generalization to Gauss-Bonnet-Maxwell-massive gravity and its stability, geometrical thermodynamics and P — V 
criticality have been investigated [0. 

One of the main problems of Maxwell’s electromagnetic field theory for a point-like charge is that there is a singular¬ 
ity at the charge position and hence, it has infinite self-energy. To overcome this problem in classical electrodynamics, 
Born and Infeld in Ref. 0 introduced a nonlinear electromagnetic field, with main motivation, to solve infinite 
self-energy problem by imposing a maximum strength of the electromagnetic field. Then, Hoffmann in Ref. (l8| 
investigated EN gravity in the presence of Born-Infeld (BI) electrodynamics. In recent two decades, exact solutions 
of gravitating black objects in the presence of BI theory have been vastly investigated 0 , 0 |. Another interesting 
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property of BI is that, BI type effective action arises in an open superstring theory and D-branes are free of physical 
singularities [21] . For a review of aspects of BI theory in the context of string theory see Ref. [|]| . Recently, there 
has been growing interest in Eddington-inspired BI gravity in context of black holes and cosmology [23} . Also, it was 
proposed that one can consider BI theory as a gravitational theory [24]. Dualization of the BI theory and some of the 
special properties of this theory have been investigated in Ref. [251 ] . 

There are several approaches for studying and obtaining critical behavior and phase transition points of black 
holes: First method is based on studying heat capacity. It was pointed out that roots and divergencies of the heat 
capacity are representing phase transition points. In other words, in place of roots and divergencies of the heat 
capacity system may go under phase transition. Another important property of the heat capacity is investigation of 
the thermal stability. Systems with positive heat capacity are denoted to be in thermally stable states. Therefore, 
the stability conditions are indicated by changes in sign of heat capacity [26|. This is known as canonical ensemble. 

In the second method, by using the renewed interpretation of cosmological constant as thermodynamical variable, 
one can modify the thermodynamical structure of the phase space [27]. One of the most important property of this 
method is the similarity of critical behavior of the black holes and ordinary thermodynamical Van der Waals liquid/gas 
systems [H]. Recently, it was pointed out that the extended phase space should be interpreted as an RG-flow in the 
space of field theories, where isotherm curves codify how the number of degrees of freedom N (or the central charge 
c) runs with the energy scale [29| . On the other hand, it was shown that variation of cosmological constant could be 
corresponded to variation of number of the colors in Yang-Mills theory residing on the boundary spacetime [30]. 

The third method is using geometrical concept for studying critical behavior. In other words, by employing a ther¬ 
modynamical potential and its corresponding extensive parameters, one can construct phase space. The divergencies 
of Ricci scalar in constructed metric are denoted as phase transition points. There are several metrics for this method 
that one can name: Weinhold [E}, Quevedo [32} and HPEM [33} which has mass as thermodynamical potential and 
Ruppeiner [dll in which entropy is considered as thermodynamical potential. These metrics are used in context of 
heat capacity. Another set of metrics was introduced in Ref. [35} which can be used in context of extended phase 
space. 

Finally, a fourth method was introduced in Ref. [35} which is based on denominator of the heat capacity. In 
this method by replacing cosmological constant with pressure in denominator of the heat capacity and solving it 
with respect to pressure, a new relation is obtained for pressure. The existence of maximum for obtained relation, 
represents the critical pressure and volume in which phase transition takes place. The behavior of system in case of 
this method is consistent with ordinary thermodynamical concepts [l5L [35} . 

The outline of the paper will be as follow. In Sec. El we introduce action and basic equations related to EN- 
Bl-massive gravity. Sec. IIIII is devoted to obtain the black hole solutions of this gravity and investigation of the 
geometrical structure of them. In the next section, we calculate conserved and thermodynamic quantities related 
to obtained solutions and check the validation of the first law of thermodynamics. In section El we study thermal 
stability of the EN-BI-massive black hole solutions in canonical ensemble. Next, we consider cosmological constant 
as pressure and study the critical behavior of the system. Then we employ the geometrical methods for investigating 
thermodynamical behavior of the system and extend this study by another method. In the last section we present 
our conclusions. 


II. BASIC EQUATIONS 

The d-dimensional action of EN-massive gravity with negative cosmological constant and a nonlinear electrody¬ 
namics is 


1 


1 
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d d Xy/—g 
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(i) 


where 77 is the scalar curvature, A = — —— ^ is the negative cosmological constant and / is a fixed symmetric 
tensor. In Eq. ©. Ci are constants and Ui are symmetric polynomials of the eigenvalues of the d x d matrix 
AC£ = yjg m f<xv which can be written as follows 

U x = [AC], U 2 = [AC ] 2 - [AC 2 ] , « 3 = [AC ] 3 - 3 [AC] [AC 2 ] + 2 [AC 3 ] , 

U A = [AC ] 4 - 6 [AC 2 ] [AC ] 2 + 8 [AC 3 ] [AC] + 3 [AC 2 ] 2 - 6 [AC 4 ] . 


Here, we want to study a particular model of nonlinear electrodynamics called BI theory which has attracted lots 
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of attentions due to its relation to effective string actions. The function L(J-) for BI theory is given as 


L{F) = 4/3 




( 2 ) 


where /3 is the BI parameter, the Maxwell invariant is F = F^ V F IJV in which F^ = d^A^ — d^A^ is the electromagnetic 
field tensor and A M is the gauge potential. 

Variation of the action m with respect to the metric tensor <? Mv and the Faraday tensor F ILV . leads to 

1 9 p \F X 

+ A — —g^ v L{F) - , M 1 + ni 2 x^u = 0, (3) 

V 1 + W 


d. 


( y^g F^ \ 


= 0, 


( 4 ) 


where G^ u is the EN tensor and x^v is the massive term with the following form 


Xnv 


-y (U,g llv - K„ v ) - y (i U 2 g V - 2+ 2 K%) - y {U 3 g^ - 3U 2 IC^ + 
6 U X K* V - 6 K%,) - y (U 4 g„ v - 4 U 3 K,^ + 12Z/ 2 AC^ - 21W,/C^ + 24K*„). 


( 5 ) 


III. BLACK HOLE SOLUTIONS IN EN-BI-MASSIVE GRAVITY 

In this section, we obtain static nonlinearly charged black holes in context of massive gravity with adS asymptotes. 
For this purpose we consider the metric of d-dimensional spacetime in the following form 

ds 2 = —f(r)dt 2 + f~ 1 (r)dr 2 + r 2 hijdxidxj, i,j = 1,2,3, ...,n , (6) 

where hi : jdxjdx : j is a (d — 2) dimension line element for an Euclidian space with constant curvature (d — 2) (d — 3 )k 
and volume Vd- 2 - We should note that the constant k, which indicates that the boundary of t = constant and 
r = constant , can be a negative (hyperbolic), zero (flat) or positive (elliptic) constant curvature hypersurface. 

We consider the ansatz metric [U| 

= diag(0,0,c 2 hij), (7) 

where c is a positive constant. Using the metric ansatz flTJ), Ui s are in the following forms [l3j 

7 . d 2 c d 2 d 3 c 2 d 2 d 3 d/ic 3 d 2 d 3 did 3 c A 

in which di = d — i. Using the gauge potential ansatz A M = h(r)5® in electromagnetic equation (j4j) and considering 
the metric (0), we obtain 


h{r) 


in which "H is the following hypergeometric function 



H = 2 F\ 


1 d 3 
2’2 ~d 2 


3 ^ 7/3 

2d 2 



( 8 ) 

(9) 


where F = and Q is an integration constant which is related to the electric charge. Also, the electromagnetic 

field tensor in d-dimensions is given by 


Ftr 


V d 2 d 3 q 
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Now, we are interested in obtaining the static black hole solutions. One may use components of Eq. m and obtain 
metric function /(r). We use the tt and x 4 x 4 components of the Eq. ©, which can be written as 


ett = {d 2 m 2 c [cir 3 + d 3 c 2 cr 2 + d 3 d 4 c 3 c 2 r + dsd^d^c^c 3 ] — 2Ar 4 


- d 2 d 3 r 2 f - d 2 r 3 f + 4/3 2 r 4 + d 2 d 3 r 2 k} - 4/3 : 


- 4/3V = 0, 


e XlXl = d 3 m 2 c [cir 3 + d 3 d 4 c 2 cr 2 + d 3 d 4 d 3 c 3 c 2 r + d 3 d<±d 3 d 3 CiC 3 ^ — 2Ar 4 

-2 d 3 r 3 f - d 3 d 4 r 2 f - r 4 f" + 4/3 2 r 4 - 0r 4 0P 2 - h' 2 + d 3 d 4 r 2 k = 0. 

We can obtain the metric function f(r), by using the Eqs. (fill) and (TH21) with the following form 


f(r) = k-^ 

f u 3 


4/3 2 - 2A 

d\d 2 




d\d 2 


dir 2d * 


2 I cci 2 , d 3 c 3 c 3 d 3 d 4 c A c 4 

+m < —r + c c 2 3-1- 5 - 


( 11 ) 


( 12 ) 


(13) 


where mo is an integration constant which is related to the total mass of the black hole. It should be noted that, 
obtained metric function (I l.'il) . satisfy all components of the Eq. m, simultaneously. 

Now, we are in a position to review the geometrical structure of this solution, briefly. We first look for the essential 
singularity(ies). The Ricci scalar and the Kretschmann scalar are 


lim R 

i ->-0 


oo, 


lim R aPlS R a ^ s 

r —>-0 


OO, 


(14) 

(15) 


and so confirm that there is a curvature singularity at r = 0. The Ricci and Kretschmann scalars are |^A and ^/tA 2 
at r —> oo. Therefore, the asymptotic behavior of these solutions are (a)dS for (A < 0 ). 

On the other hand, in the absence of massive parameter (m = 0), the solution (1131) reduces to an d-dimensional 
asymptotically adS topological black hole with a negative, zero or positive constant curvature hypersurface in the 
following form 


f(r) 



aitt 2 \d1d2 


1^2, 4 d 2 q 2, H 

I 2 ) dir 2d3 


(16) 


In order to study the effects of the EN-BI-massive gravity on metric function, we have plotted various diagrams 

(Figs. m-m). 

By considering specific values for the parameters, metric function has different behaviors. Depending on the choices 
of the parameters, EN-BI-massive black holes can behave like Reissner-Nordstrom black holes. In other words, these 
black holes may have two horizons, one extreme horizon and without horizon (naked singularity) (see Fig. (T] for more 
details). On the other hand, by adjusting some of the parameters of EN-BI-massive black holes, we encounter with 
interesting behaviors. The solutions may have three or higher horizons (Figs. [2] and [3]). The existence of three or 
higher horizons for black holes is due to the presence of massive gravity [lji, [l(| . In addition to the significant effects 
of massive term, we should note that the nonlinearity parameter can affect on the number of horizons. In addition, /3 
can change the type of singularity. In other words, depending on the parameters, one can find a /3 C in which singularity 
is spacelike for /3 < /3 C , and it would be timelike for /3 > /3 C (see [2(J for more details). 

Now, we give a brief discussion regarding Carter-Penrose diagrams. In order to study the conformal structure of 
the solutions, one may use the conformal compactification method through plotting the Carter-Penrose (conformal) 
diagrams. As we mentioned before, depending on the value of nonlinearity parameter, /3, one may encounter with 
timelike or spacelike singularity. Penrose diagrams regarding to timelike singularity was discussed in [161 ] (see conformal 
diagrams in 0). Here we focus on special case in which singularity is spacelike (/3 < j3 c ). In other words, the 
singularity of this nonlinearly charged black holes behaves like uncharged Schwarzschild solutions (see Fig. m- This 
means that, although massive and nonlinearity parts of the metric function can change the type of singularity and 
horizon structure of black holes, they does not affect asymptotical behavior of the solutions. Drawing the Carter- 
Penrose diagrams, we find the causal structure of the solutions are asymptotically well behaved. 
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FIG. 1: f(r) versus r for A = —0.8, q = 1, 8 = 0.6, m = 1.4, c = —0.8, ci = 2, C3 = —4, C4 = 0, k = 1, and d = 4. 

Left diagram for m 0 = 5, C2 = 1.00 (dashed line), C2 = 1.22 (continues line) and C2 = 1.35 (dotted line). 

Right diagram for C2 = 1.40, mo = 5.80 (dashed line), mo = 5.64 (continues line) and mo = 5.5 (dotted line). 
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FIG. 2: f(r) versus r for A = —1, q = 0.5, ft = 7, m = 0.5, c = 0.4, ci = —40, C2 = 60, C 3 = 1, C4 = 0, k = 1 and d = 6. 
diagrams for mo = 1.75 (dashed line), mo = 1.67 (continues line) and mo = 1.61 (dotted line). 

IV. THERMODYNAMICS 

In this section, we calculate the conserved and thermodynamic quantities of the static EN-BI-massive black hole 
solutions in d-dimensions and then check the first law of thermodynamics. 

By using the definition of Hawking temperature which is related to the definition of surface gravity on the outer 
horizon r + , one can find 

T = -— T [ciri + d^cr 2 , + dzdAC^c 2 r+ + d3d4dsC4C 3 l + - —— -h - — y/l + T+ 

47T r+ L J 27ra2 47rr + nd2 




(17) 
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FIG. 3: Metric functions and Carter Penrose diagrams for the asymptotically adS black holes with spacelike singularity. Three 
horizons (continuous line of metric function and related Carter Penrose diagram in right-up panel), two horizons which inner 
one is extreme (dotted line of metric function and related Carter-Penrose diagram in left-down panel) and two horizons which 
outer one is extreme (dashed line of metric function and related Carter-Penrose diagram in right-down panel). 


where T + = d2d lj 0 . The electric charge, Q , can be found by calculating the flux of the electric field at infinity, yielding 
f3 2 r , 2 


Vd 2 Vd2<h 
Q = —^ q - 


(18) 


In order to obtain the entropy of the black holes, one can employ the area law of the black holes. It is a matter of 
calculation to show that entropy has the following form [36] 


c _ Y^± r d 2 

* “ ~ r + ’ 


(19) 


It was shown that by using the Hamiltonian approach or counterterm method, one can find the mass M of the 
black hole for massive gravity as 


„ T d 2 V d2 

M = “T7— m o> 

167T 

in which by evaluating metric function on the horizon (/ (r = r+) = 0), one can obtain 
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where H+ = 2^1 


1 _ds_ 
2’ 2 d 2 


3d- 


7/3 


2d 2 


,-r. 


-)■ 


It is notable that, U is the electric potential, which is defined in the following form 


U = I™ - |r-r + = \H L \ n+. (22) 

V ^3 tv* 

Having conserved and thermodynamic quantities at hand, we are in a position to check the first law of thermody¬ 
namics for our solutions. It is easy to show that by using thermodynamic quantities such as charge (usd, entropy dD 
and mass urn with the first law of black hole thermodynamics 


dM = TdS + UdQ, 


(23) 


we define the intensive parameters conjugate to S and Q. These quantities are the temperature and the electric 
potential 


T = 



and 



(24) 


which are the same as those calculated for temperature m and electric potential (l22l) . 


V. HEAT CAPACITY AND STABILITY IN CANONICAL ENSEMBLE 

Here, we study the stability conditions and the effects of different factors on these conditions. The stability 
conditions in canonical ensemble are based on the signature of the heat capacity. The negativity of heat capacity 
represents unstable solutions which may lead to following results: unstable solutions may go under phase transition 
and acquire stable states. This phase transition could happen whether when heat capacity meets a root(s) or has 
a divergency. Therefore, the roots of regular numerator and denominator of the heat capacity are phase transition 
points. In the other scenario, the heat capacity is always negative. This is known as non-physical case. But there is 
a stronger condition which is originated from the temperature. The positivity of the temperature represents physical 
solutions whereas its negativity is denoted as non-physical one. Therefore, in order to getting better picture and 
enriching the results of our study, we investigate both temperature and heat capacity, simultaneously. 

The heat capacity is given by 
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( d 2 AI \ 


(dT\ ■ 
V dS ) Q 


Considering Eqs. C3 and (fTTJl) . it is a matter of calculation to show that 
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Trc^rJ 
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7rd2r , + 

4 d 3 q 2 (d\ + T 
Trd 2 r^ 7 ^/l + 


[3d 4 d 5 c 4 c 2 + 2 d 4 C 3 cr + + C 2 Cr^] 


(26) 


In order to study the effects of different parameters on stability conditions and temperature, we have plotted various 
diagrams (Figs. |4]-[8|). 

Interestingly, in the absence of the massive parameter (Fig. [4] right panel), temperature starts from —00 and it is 
only an increasing function of horizon radius with a root. Therefore, we have two regions of physical and non-physical 
solutions. Adding massive gravity could modify the behavior of the temperature into an U shape diagram starting 
from +00 without any root. The extremum is an increasing function of massive parameter (Fig. [4] right panel) and 
dimensions (Fig. [5] right panel), whereas, it is a decreasing function of k (Fig. [6] right panel) and electric charge 
(Fig. [7] right panel). The only exception for this behavior is for strong nonlinearity parameter. Interestingly, for 
large values of nonlinearity parameter, a massive-less like behavior is observed (Fig. [5] right panel). In other words, 
the temperature starts from —00 but the effect of the massive could be seen through two extrema. The root of 
temperature and smaller extremum are increasing functions of /3 and larger extremum is a decreasing function of it. 
Another interesting property of temperature is the effect of dimensions. Studying Fig. [8] (right panel) shows that the 
temperature for every two sets of dimensions will coincide with each other. In other words, there are places in which 

















FIG. 4: For different scales: Cq (left panel) and T (right panel) versus r+ for q = 1, A = —1 c = ci = C 2 = C 3 = 2, a = 0, 
f3 = 0.5, d = 5 and k = 1; m = 0 (continues line), m = 0.25 (dotted line), m = 0.35 (dashed line) and m = 0.40 (dashed-dotted 
line). 



FIG. 5: For different scales: Cq (left and middle panels) and T (right panel) versus r+ for q = 1, A = —1, c = c\ = C2 = C3 = 2, 
C4 = 0, m = 0.30, d = 5 and k = 1; /3 = 2 (continues line), /3 = 3 (dotted line), /3 = 4 (dashed line) and (3 = 5 (dashed-dotted 
line). 


despite differences in dimensions, two black holes with two different dimensions and same values for other parameters 
will have same temperature in a special r+. 

In absence of massive gravity, black holes could acquire temperature from zero to + 00 , whereas adding massive, 
will cause the black holes never acquire some temperature. This effect is vanished in case of large nonlinearity 
parameter. In other words, the strength of nonlinearity parameter has opposing effects to massive’s ones. Also, the 
U shape diagram indicates that for every temperature that black holes can acquire two horizons exist except for the 
extremum. Therefore, considering Hawking radiation, one is not able to recognize the size of these black holes by 
measuring their Hawking radiation. It is worthwhile to mention that extrema and root(s) of temperature are phase 
transition points of heat capacity. 

Regarding stability, it is evident that in absence of the massive gravity, there exists a region of the instability 
which is located where the temperature is negative. Therefore, this is a non-physical solution (Fig. [T] left panel). 
Interestingly, by adding massive gravity, the non-physical region is vanished and heat capacity acquires divergence 
point without any root. Before divergence point, the heat capacity is negative. Therefore, in this region black holes 
are unstable. In divergence point, black holes go under phase transition of smaller unstable to larger stable black 
holes. The divergence point is an increasing function of massive parameter (Fig. [3] left panel) and dimensions (Fig. [ 8 ] 
left panel), whereas, it is a decreasing function of k (Fig. [ 6 ] left panel) and electric charge (Fig. [7] left panel). 
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FIG. 6: For different scales: Cq (left and middle panels) and T (right panel) versus r+ for q = 1, A = —1, c = ci = C2 = C3 = 2, 
C4 = 0, m = 0.4, d = 5 and /3 = 0.5; k = 1 (continues line), fc = 0 (dotted line) and k = — 1 (dashed line). 




FIG. 7: For different scales: Cq (left panel) and T (right panel) versus r+ for A = —1, c = ci = C 2 = C 3 = 2, ci = 0, /3 = 0.5, 
m = 0.4, d = 5 and k = 1; q = 0 (continues line), q = 0.5 (dotted line), q = 1 (dashed line) and q = 2 (dashed-dotted line). 


Interestingly, in strong nonlinearity parameter, the mentioned behavior is modified. In this case black holes enjoy 
one root and two divergence points. Before root and between two divergencies, heat capacity is negative and between 
root and smaller divergence point and after larger divergence point, heat capacity is positive. According to thermo¬ 
dynamical concept, systems go under phase transition to acquire stable states. Therefore, following phase transitions 
take place: non-physical unstable to physical stable (in root), large unstable to smaller stable (in smaller divergence 
point) and smaller unstable to larger stable black holes (in larger divergency). Root and smaller divergence point 
are increasing functions of /3 (Fig. [5] left panel), whereas, larger divergency is a decreasing function of it (Fig. [5] 
middle panel). It is worthwhile to mention that larger divergency is not highly sensitive to variation of nonlinearity 
parameter. 

Comparing obtained results for heat capacity (regarding phase transitions) and the behavior of the temperature, 
one can see that larger to smaller phase transition takes place at maximum (compare Fig. [G] left diagram with right) 
and smaller to larger one happens at minimum (compare Fig. ^middle diagram with right) of temperature. Therefore, 
one is able to recognize the type and number of phase transition by only studying temperature’s diagrams. 
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FIG. 8: For different scales: Cq (left and middle panels) and T (right panel) versus r+ for q = 1, A = —1, c = ci = C2 = C3 = 2, 
C4 = 0, P = 0.5, m = 0.4 and k = 1; d = 5 (continues line), d = 6 (dotted line) and d = 7 (dashed line). 


VI. P - V CRITICALITY OF CHARGED BLACK HOLES IN EN-BI-MASSIVE GRAVITY 


Now, we are in a position to study the critical behavior of the system through phase diagrams. Using the renewed 
interpretation of the cosmological constant as thermodynamical pressure, one can use following relation to rewrite 
thermodynamical relations of the solutions in spherical horizon [281 ] 


P= — 


A 

8w ’ 


which results into following conjugating thermodynamical variable corresponding to pressure 

v=(-) ■ 

\dPJs,Q 


(27) 


(28) 


Due to existence of the pressure in obtained relation for total mass of the black holes, one can interpret the total 
mass as thermodynamical quantity known as Enthalpy. This interpretation will lead to the following relation for 
Gibbs free energy [28[ 


G = H-TS = M-TS. 

Now by using Eqs. (EU1) and (1271) with the relations of volume and Gibbs free energy (Eqs. 


and 


rr d\ 

F= v+- 


(29) 
), one finds 

(30) 


and 


G = 


-P 


G?lC?2 


m 2 c 2 r([ 5 

167T 


(3d 3 d4C4C 2 + 2d 3 C3cr + + C2r\) 


dh 2 u+ 

27rdiri 3 


2„d 

+ yi+rv+ ' + 


47T<ii<i2 


167T 


(31) 


Obtained relation for volume indicates that volume of the black holes is only a function of the topology of the 
solutions and independent of electrodynamics and gravitational extensions, directly. 

In order to obtain critical values, one can use P—V diagrams. In other words, by studying inflection point properties 
one can obtain critical values in which phase transitions may take place. Therefore, one can use 


dP 

dr + 


- ( d2p \ 


dr\ Jr 


= 0. 


(32) 
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FIG. 9: P — r+ (left), T — r+ (middle) and G — T (right) diagrams for /? = 0.5, q = 1, m = 0.1, c = ci = C 2 = C 3 = 0.2, ca = 0 
and d = 5. 

P — r+ diagram, from up to bottom T = 1.17k, T = T c and T = 0.971, respectively. 

T — r+ diagram, from up to bottom P = 1.17k, P = P c and P = 0.9 P c , respectively. 

G — T diagram for P = 0.5P C (continuous line), P = P c (dotted line) and P = 1.5P C (dashed line). 


Considering obtained values for temperature 03 and pressure m, one can obtain pressure as 

P = - Y [ c i r + + rf 3C2 cr\ + d 3 d 4 c 3 c 2 r+ + d 3 d 4 d 5 c 4 c 3 ] - (V 1 + r + ~ l) • (33) 

Now, by considering Eq. ( 15 ^ 1 ) with obtained relation for pressure ( 1551 ) . one can obtain two relations for finding critical 
quantities. Due to economical reasons, we will not present them. Regarding the contribution of electromagnetic part, 
it is not possible to obtain critical horizon analytically, and therefore, we use numerical method. Considering the 
variation of /3 and massive parameter, one can draw following tables 


m 

r c T c P c 

0.000000 

1.8264628 0.1334354 0.02200146 0.3011558 

0.100000 

1.8263848 0.1334835 0.02200495 0.3010822 

1.000000 

1.7953522 0.1382092 0.02233685 0.2901581 

5.000000 

1.6278897 0.2562956 0.03187871 0.2024811 

10.000000 0.7052643 0.7329062 0.13205092 0.1270705 

Table (1): q = 

1, P = 0.5, ci = C 2 = c 3 = 0.2, c 4 = 0 and d = 


P r c _ Tc _ Pc 

1.000000 1.7819632 0.1693410 0.0296043 0.3115245 
2.100000 1.8174387 0.1677996 0.0290170 0.3142835 
3.500000 1.8235114 0.1675323 0.0289161 0.3147386 
4.000000 1.8256045 0.1674399 0.0288813 0.3148944 
5.000000 1.8265678 0.1673974 0.0288653 0.3149659 
Table (2): q = 1, m = 0.1, Ci = C 2 = c 3 = 2, c 4 = 0 and d = 5. 


In addition, we plot following diagrams (Figs. [9] - fl2l) to investigate that obtained values are the ones in which 
phase transition takes place or not. 

The formation of swallow tail in G — T diagrams for pressure smaller than critical pressure (Figs. [5] and |TT] right 
panels), subcritical isobars in T — r+ diagrams for critical pressure (Figs, [ill and [Til middle panel) and isothermal 
diagrams in case of critical temperature in P — r+ diagrams (Figs. [9] and [lT] left panels), show that obtained values 
are critical ones in which phase transition takes place. 
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FIG. 10: P — r+ (left), T — r+ (middle) and G — T (right) diagrams for /3 = 0.5, q = 1, c = ci = C 2 = C 3 = 0.2, C 4 = 0, d = 5, 
m = 0 (continuous line), m = 1 (dotted line) and to = 5 (dashed line). 

P — r+ diagram for T = T c , T — r+ diagram for P = P c and G — T diagram for P = 0.5 P c . 



FIG. 11: P — r+ (left), T — r+ (middle) and G — T (right) diagrams for /3 = 2, q = 1, m = 0.1, c = ci = C 2 = C 3 = 2, C 4 = 0 
and d = 5. 

P — r+ diagram, from up to bottom T = 1.1T C , T = T C and T = 0.9T C , respectively. 

T — r+ diagram, from up to bottom P = 1.1P C , P = P c and P = 0.9P C , respectively. 

G — T diagram for P = 0.5P C (continuous line), P = P c (dotted line) and P = 1.5P C (dashed line). 


0.0310 


0.0305 


0.0300 


0.0295 


0.0290 


0.0285 





FIG. 12: P — r+ (left), T — r+ (middle) and G — T (right) diagrams for m = 0.1, q = 1, c = ci = C 2 = C 3 = 0.2, C 4 = 0, d = 5, 
(H — 1 (continuous line), /3 = 2 (dotted line) and /? = 3 (dashed line). 

P — r+ diagram for T = T c , T — r+ diagram for P = P c and G — T diagram for P = 0.5P C . 
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It is evident that critical pressure (Fig. [TO] left panel) and temperature (Fig. [TO] middle panel) are increasing 
functions of the massive parameter, whereas the critical horizon (Fig. [TO] left and middle panels) and universal ratio 
of are decreasing functions of this parameter. 

It is worthwhile to mention that length of subcritical isobars (which is known as phase transition region) is a 
decreasing function of massive parameter (Fig. [TO] middle panel). In opposite, the size of swallow tail and the energy 
of different phases are increasing functions of m (Fig. [TO] right panel). 

Interestingly, the effects of variation of nonlinearity parameter is opposite of massive parameter. In other words, 
critical pressure (Fig. [ 12 ] left panel), temperature (Fig. [12] middle panel) and the size of swallow tail (Fig. IT 2 l left 
panel) are decreasing functions of /?, whereas, the critical horizon radius (Fig. [12] left and middle panels), length 
of subcritical isobars (Fig. [12] middle panel) and universal ration of are increasing functions of nonlinearity 
parameter. 

It should be pointed that the length of subcritical isobars affects single regions of different states which in our cases 
are smaller and larger black holes. In other words, increasing the length of subcritical isobars (phase transition region) 
decreases the single state regions. 


A. Neutral Massive black holes 


In this section, by cancelling the electric charge (q = 0), we will study the critical behavior of the system. Previously, 
it was shown that Schwarzschild black holes does not have any phase transition in context of extended phase space. 
Now, we are investigating the effects of massive gravity in case of EN-massive gravity. Using obtained relation for 
calculating critical horizon radius in previous part and setting q = 0 , one can find following relation for calculating 
critical horizon radius 


m 2 ( Qd^d^CiC 4 + 3 c? 4 C 3 C 3 r + + c 2 c 2 r 2 [_) + r\ = 0. (34) 

It is a matter of calculation to show that this relation has following roots which are critical horizon radii 

me 2 ( 3 d 4 ?nC 3 C± \f~3di [ 8 ds (1 + C 2 C 2 to 2 ) C 4 — 3cLito 2 c§c 2 ] 1 
r c = -^- oTTH - 2 - 2 T-- ■ (35) 

2(1+ C 2 C“TO z ) 

Obtained relation shows that in absence of massive gravity, critical horizon radius will be zero which is not of our 
interest. This result consistent with Schwarzschild case. Now, for the simplicity, we consider the case of C 4 = 0. This 
leads into following critical horizon radius 

Sd^mcsc 


T nr. — 


1 + C 2 C 2 TO 2 


(36) 


It is evident that for the cases of d = 4 and d > 4 with vanishing C 3 , the critical horizon radius will be zero. 
Therefore, there is no phase transition for these black holes. Interestingly for case of d > 4, the condition for having 
a positive critical horizon radius will be C 3 < 0 and 1 + C 2 C 2 to 2 > 0. By employing obtained value for critical horizon 
radius, one can find critical temperature and pressure in the following forms 


( 3^40403 — d 3 c|) m 4 c 4 — d% ( 2c 2 m 2 c 2 + l) 
127rd4?7? 2 C3C 3 


(37) 


_ (c 2 to 2 c 2 + I ) 3 d- 2 d 3 
cc 4327T d 2 m 4 c|c 6 

Considering obtained values, one can show that following equality is hold 

Pc C r C c _ _ (c 2 to 2 c 2 + if d 2 d 3 _ 

T cc 12 ( 3 d 4 CiC 3 — d 3 c 2 ) m 4 c 4 — 12^3 ( 2 c 2 m 2 c 2 + 1) ’ 


(38) 


(39) 


which shows that in this case, Pc £ rac is a function of massive parameter and coefficients. Using obtained critical values 
(Eqs. (1361) . (I3T1) and ([38]) ) with Eqs. lHTjl . ([27jl . (l29l) . (l33l) and setting q = 0, we plot following diagrams for 5 and 6 
dimensions (Figs. [151 and [TH). 

In Ref. [37|, it was shown that in context of neutral Gauss-Bonnet black holes, no phase transition is observed 
in 6 -dimensions. Here, the extension of the massive gravity enables the black holes to enjoy the existence of second 
order phase transition in 6 -dimensions. Also, we observed that contrary to Gauss-Bonnet case, p ^ rcc is a function of 
massive gravity. 
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FIG. 13: P — r+ (left) and G — T (right) diagrams for m = 5, c = 0.2, ci = C2 = 2, C3 = —2, a — 0 and d = 5. 
P — r+ diagram, from up to bottom T = 1.1T C , T = T c and T = 0.9 T c , respectively. 

G — T diagram for P = 0.5P C (continuous line), P = P c (dotted line) and P = 1.5 P c (dashed line). 



FIG. 14: P — r+ (left) and G — T (right) diagrams for m = 5, c = 0.2, ci = C 2 = 2, C 3 = —2, a — 0 and d = 6. 
P — r+ diagram, from up to bottom T = 1.1T C , T = T c and T = 0.971, respectively. 

G — T diagram for P = 0.5P C (continuous line), P = P c (dotted line) and P = 1.5P C (dashed line). 


VII. GEOMETRICAL PHASE TRANSITION IN CONTEXT OF HEAT CAPACITY AND EXTENDED 

PHASE SPACE 


In this section, we employ the geometrical concept for studying thermodynamical behavior of the obtained solutions. 
In order to do so, we employ HPEM method. In this method, the thermodynamical phase space is constructed by 
considering mass of the black holes as thermodynamical potential. By doing so, the components of the phase space 
will be extensive parameters such as electric charge, entropy and etc. The general form of HPEM metric is [33j] 


ds 


2 

New 


SM s 


(rrn d 2 M \ 

1 iJ -i=2 dx'f ) 



(40) 


where Mg = dM/dS , Mgs = d 2 M/dS 2 and \i (Xi 7^ S) are extensive parameters which are components of phase 
space. Now, we will investigate whether the phase transition points that were obtained in section (IV) coincide with 
all divergencies of the Ricci scalar of HPEM metric. For economical reasons, we only plot diagrams correspond to 
variation of massive and nonlinearity parameters. To do so, we use Eqs. GHD, (El) and (1201) with HPEM metric (Eq. 
l40l) . This leads into following diagrams (Fig. fl5l) . 

It is evident that employed metric has consistent results with what were found in case of heat capacity (compare 
Figs. HI and l5l with fl5l) . In other words, the divergencies of the Ricci scalar are matched with phase transition points 
of the heat capacity. An interesting characteristic behavior of the diagrams is the different divergencies for different 
types of phase transition. In case of larger to smaller black holes phase transition, the divergency of the Ricci scalar 
is toward +00 (compare Fig. [5] left panel with Fig. E] middle panel), whereas in case of smaller to larger phase 
transition, the divergency is toward —00 (compare Fig. [5] middle panel with Fig. 1151 right panel). This specific 
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FIG. 15: For different scales: TZ versus r+ diagrams for q = 1 , c ■* ci = C2 » C3 = 2, C4 = 0, d = 5 and k = 1 . 

left: (3 = 0.5, m = 0 (continues line), m = 0.25 (dotted line), m = 0.35 (dashed line) and m = 0.40 (dashed-dotted line). 

middle and right: m = 0.3, (3 = 2 (continues line), (3 = 3 (dotted line) and (3 = 4 (dashed line). 



FIG. 16: For different scales: TZ (continuous line), Cq (dashed line) diagrams for q = 1, c = ci = ci = C3 = 2, C4 = 0, (3 = 0.5, 
d = 5, k = 1 and m = 0.1. 

P = 0.9P C left and middle panels, P = P c left and right panels and P = 1.1P C left panel. 


behavior enables us to recognize the type of phase transition independent of heat capacity. 

Next, we employ another geometrical metric for studying the critical behavior of the system in context of extended 
phase space. In this metric, Due to consideration of the cosmological constant as thermodynamical pressure, we 
have three extensive parameters; electric charge, entropy and pressure. In order to construct phase space we employ 
following metric (35j 


ds 2 = S^-(-M S sdS 2 + M QQ dQ 2 + dP 2 ) . (41) 

Considering Eqs. (ED, ED, ED and ED, we plot following diagram (Fig. flfill with respect to Fig. [9] 

Due to existence of a root for heat capacity, in all plotted diagrams, a divergency is observed (Fig. E]left panel). 
It is evident that for pressures smaller than critical pressure, system goes under two phase transitions with different 
horizon radii (Fig. [16] middle panel). This is consistent with what was observed in studying T—r+ diagrams of Fig. [9] 
On the other hand, for critical pressure system goes under a phase transition. The place of this divergency is exactly 
located at the critical horizon which is obtainable through T — r_|_ diagrams of Fig. [9] (Fig. [TB] right panel). Finally, 
for pressures larger than critical pressure no phase transition is observed and the behavior of Ricci scalar will be what 
is plotted in Fig. [16] (left panel). These results are consistent with ordinary thermodynamical concepts and indicates 
that these three pictures (phase diagrams, heat capacity and geometrical thermodynamics) are in agreement. 
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FIG. 17: P versus r + diagrams for q = 1, c = Ci = C 2 = C 3 = 2, C 4 = 0, d = 5 and k = 1. 

left panel:/! = 0.5 and m = 0 (bold-continues line), P = 0.02200146 (continues line), m = 5 (bold-dashed line) and P = 
0.03187871 (dotted line). 

right panel: m = 0.1 and [3 = 1 (bold-continues line), P = 0.0296043 (continues line), [3 = 5 (bold-dashed line) and P = 
0.0288653 (dashed line). 


VIII. HEAT CAPACITY AND CRITICAL VALUES IN THE EXTENDED PHASE SPACE 

The final section of this paper is devoted to calculation of the critical pressure in extended phase space by using 
denominator of the heat capacity. It was shown that one can calculate critical pressures that were obtained in section 
( V) by using denominator of the heat capacity [35}. To do so, one should replace the cosmological constant in 
denominator of the heat capacity m with its corresponding pressure m ■ Then, solve the denominator of the heat 
capacity with respect to pressure. This will lead into following relation 

P = 2 ^^4 — (3d 4 d 5 c 4 c 2 + 2d 4 c 3 cr + + c 2 r \) 

(\/l + r+ — l) P 2 d 2 d^ 

47Ta/ 1 + T+ 167rr+ 

Obtained relation for pressure is different from what was obtained through use of temperature (1351) . In this relation, 
the maximum(s) of pressure and its corresponding horizon radius are critical pressure and horizon radius in which 
phase transition takes place. Now, by using indicated values in table 1 and Eq. HUD, we plot following diagram (Fig. 

[ 13 - 

It is evident that obtained maximums are critical pressures in which phase transitions take place. The thermody¬ 
namical concept that was mentioned in last section (pressure being smaller than critical pressure leads to existence 
of two phase transitions and for pressures larger than critical pressure no divergency is observed) is also hold in case 
of this approach. In other words, this approach is an additional method for studying critical behavior of the system 
and the results of this approach is consistent with GTs, heat capacity and extended phase space. 


d 2 d 2 q 2 


8 tt r 2 + d2 




(42) 


IX. CONCLUSIONS 

In this paper, we have considered EN-massive gravity in presence of BI nonlinear electromagnetic field. It was 
shown that considering this configuration leads to modification of the number and place of horizons that black holes 
can acquire. In other words, cases of multiple horizons were observed with different phenomenologies. Next, conserved 
and thermodynamical quantities were obtained and it was shown that first law of thermodynamics hold for these black 
holes. 

Next, we studied the thermodynamical behavior of the system. It was shown that temperature in the presence and 
absence of massive gravity presents different behaviors. Adding massive put limitations on values that temperature can 
acquire, while, there was no limitation for temperature in the absence of it. Interestingly, this behavior was modified 
in the presence of strong nonlinearity parameter. In strong nonlinearity parameter, the behavior of temperature 
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returned to a massive-less like behavior but the effects of the massive were observed in existences of extrema. It was 
also seen, that in case of different dimensions, for each pair of dimensions, one can find a point in which temperature 
for both dimensions are equal. 

Regarding the stability, it was seen that in the presence of massive gravity, black holes enjoy one phase transition of 
the smaller unstable to larger stable. The phase transition was related to the divergency of heat capacity. Then again, 
in strong nonlinearity parameter, this behavior was modified. In this case, black holes had three phase transitions of 
smaller non-physical unstable to larger physical stable (in place of root), larger unstable to smaller stable (in place of 
smaller divergency) and smaller unstable to larger stable (in place of larger divergency). 

Clearly, one can conclude that nonlinear electromagnetic field has an opposing effect comparing to massive gravity. 
Strong nonlinearity parameter modifies the effects of the massive gravity and return the system to the massive-less 
like behavior, although the effects of massive still observed through extrema. 

It was pointed out that at maximums of the temperature, larger unstable to smaller stable and at minimums, smaller 
unstable to larger stable phase transitions take place. Therefore, studying temperature provides an independent 
picture for studying phase transitions and stability of the solutions. 

Next, we extended phase space by considering cosmological constant as thermodynamical variable known as pressure. 
It was shown that volume of the black holes is independent of generalization of the electromagnetic field and extension 
of the massive gravity. Obtained values where critical points in which phase transitions took place. It was shown that 
the effects of variation of nonlinearity parameter was opposite of the massive parameter. In other words, these two 
factors put restrictions on each others effects. 

Interestingly, in Ref. jlH ] variation of massive gravity highly modified the critical temperature and pressure. In 
case of obtained solutions in this paper, the modification was not as considerable as what was observed in case of 
Gauss-Bonnet-Maxwell-massive black holes. This shows that generalization of electromagnetic field puts stronger 
restrictions on the effects of the massive gravity. In other words, in order to have stronger control over contributions 
of the massive gravity one should increase the nonlinearity of the electromagnetic sector. 

In addition, a study in context of neutral solutions was conducted. It was shown that due to contribution of the 
massive gravity, the chargeless solutions of this gravity also enjoy the existence of phase transition. In other words, 
black holes in EN-massive gravity go under phase transitions in extended phase space. Also, it was shown that ratio of 
was a function of massive gravity. It was also pointed out that in 6-dinrensions, contrary to case of Gauss-Bonnet 
black holes, these black holes enjoy second order phase transition. In addition, it was shown that in case of d = 4, no 
phase transition for massive black holes is observed. 

Next, geometrical approach was used for studying critical behavior of the system in context of heat capacity and 
extended phase space. It was shown that employed metrics for both cases have consistent results and follow the 
concepts of ordinary thermodynamics. The characteristic behavior of divergencies in Ricci scalar of the geometrical 
thermodynamical metrics, enabled us to recognize the type of phase transition (smaller to larger or larger to smaller). 

Finally, another method which was based on denominator of the heat capacity was used to calculate critical pressure 
and horizon radius. It was shown that this method has consistent results with extended phase space and follow the 
concepts of ordinary thermodynamics. In other words, this method provides an independent approach for investigating 
critical behavior of the system. 

Due to generalization of Born-Infeld for electromagnetic sector of the solutions, it will be worthwhile to study 
the effects of this generalization on conductivity of these black holes and their corresponding superconductors phase 
transition. Specially, it will be interesting to see how this generalization will affect the interpretation of graviton as 
lattice and the Drude like behavior. Also, it will be worthwhile to study the metal-insulator transition in context of 
these solutions. 
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